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1 . Introduction 

A triangular  matrix  reveals  its  eigenvalues  on  the  main  diagonal.  By 
Schur's  lemma  any  square  matrix  is  unitarily  similar  to  an  upper  triangular 
matrix  with  the  eigenvalues  arranged  in  any  desired  order  alonq  the  diaqonal. 
In  practice  the  QR  algorithm  in  real  arithmetic  produces  a block  triangular 
matrix  in  which  the  eigenvalues  are  likely  to  be  in  monotone  decreasing 
order  by  absolute  value  down  the  diagonal.  However  this  monotonic ity 
cannot  be  guaranteed  and  for  some  purposes  the  ordering  bv  absolute  value 
is  not  what  is  wanted. 

The  problem  which  we  address  here  is  to  find  some  simple  orthogonal 
similarity  transformations  which  have  the  effect  of  exchanging  two  diagonal 
elements  (or  blocks)  while  preserving  block  triangular  form.  Actually  we 
will  show  only  how  to  swap  adjacent  blocks  and  so  the  exchange  of  distant 
blocks  must  be  accomplished  by  a succession  of  adjacent  swaps, 

Although  the  cost  of  such  a swap  is  snail  it  is  not  negl igibleT^Kf  an 
n x n matrix  (p+q^n  multiplications  are  needed  to  swap  adjacent  diagonal' 
blocks  of  orders  p and  q. 


2.  Ruhe’s  Trick 


For  any  real  0 and  s * sin  O,  c * cos  0 
-s  c 

is  an  orthogonal  matrix  representing  a 
Observe  that 


the  symmetric  matrix 
reflection  of  the  plane. 


-s 

c 

ct.  e 

-s  c 

c 

' a^^-Bsc+a^c  , -o^sc-Bs  -M^sc 
,2  2 , , 2 

c 

s 

0 

c s 

-a.sc+Bc  +a^sc  , aijC  +PSC+U2S 

The  new  matrix  is  upper  triangular  if  and  only  if 


1 


c ( Be  - (o^-Os)  3 0 . 


4 

The  choice  c * 0 represents  no  change,  the  choice 

tan  0 * s/c  * ty(ou-ap) 

results  In  an  exchange  of  a,  and  o?.  The  new  (1,?)  element  Is 

2 

-s[BS  + (a1-a2)c]  * -s[fls*tfc  /s|  - -P  . 

Now  suppose  that  is  the  (j,j)  element  of  an  n " n upper  triangular 
matrix.  The  plane  reflection  Indicated  above,  effected  in  the  (%i,j*ll 
coordinate  plane,  will  swap  and  o0.  Postmultiplication  affects 
columns  j and  j+1  while  premultlplication  affects  rows  j and  ,1+1. 
This  requires  4 ( n-2 ) multiplications.  To  keep  the  angle  0 in 
(-n/2,n/2)  we  define 

d - . 

c * lo.-Opl/d 
s » s1gn(n^-«^)/d  . 

Note  that  when  H a D the  transformation  merely  exchanges  the  two  rows  and 
the  corresponding  pair  of  columns. 

3 . The  General  Case 
Consider  the  reduced  matrix 

A1  8 A1  Is  pxp  , 

,0  A2  Is  q * q . 

We  seek  an  orthogonal  similarity  transformation  which  swaps  A1  and  A„. 

In  general  this  Is  not  possible;  fortunately  we  can  achieve  a form  which 
Is  as  useful  as  exchanging  A,  and  A0.  We  denote  b\  7 1 the  transpose 


. ZST’r  . rsST  » Tm 


of  any  matrix  2.  A partitioned  matrix 


_ST  c 

l2 


C,  Is  p x p , 


l C1  s2  J C2  is  qxq  , 

Is  orthogonal  If,  and  only  If,  the  following  relations  hold: 

«>  C1C>S2S2-  'p-  VXC,  • 

/ o \ cTc  . r _ t - . c^c 


(1) 

^1^1  + ^2^2  s 

'p* 

(2) 

S1S1+C2C2  = 

V 

(3) 

- C]S1  +S2C2  * 

°p.q 

(4) 

~ C2 + c^s2  = 

°q.p 

Note  that  if  c|  = = C2  then  we  can  take  * S2,  however  this 

is  not  always  advantageous. 

We  seek  an  orthooonal  matrix  of  the  form  shown  above  such  that 


-S]  C2  A]  B 


c1  s2  j i 0 a2 


a2  b f -s;  c2 


0 A,  Jl  C,  S2 


On  equating  the  (2,1)  and  (2,2)  blocks  on  each  side  of  the  equation  we  find 

(5)  C1A1  * A]C1  (also  A2C2  = Cpy  , 

(6)  C.B  + S-A,  * A,S,  . 


C]B + S2A2  = A^S2  . 


When  C,  is  Invertible  (more  on  this  below)  then  (6)  can  be  rewritten  as 

8 + C1  S2A2  = C1  k^S2 


A1C;'S2  , by  (5)  . 


z 


6 


r 


We  now  let  the  p * q matrix  Cj^S2  ■ X/f.,  where  f.  is  a positive  constant 
at  our  disposal,  and  substitute  into  the  equation  above  to  get 

(7)  A^X  - XA?  - f.B  . 

In  order  to  obtain  from  X we  pre-  and  post-multiply  the  first  ortho- 
gonality relation  (1)  appropriately  and  Invert  to  find 

Ip  + XX '/t  a c1  'c^ 
or 

(8)  (C1/C)T(C1/C)  a (l/^XX7)-1  W,  . 

Using  (3)  we  find  that  X/C.  also  equals  and  by  usinq  (?)  we  obtain 

(9)  (c2/r.)T(c2/0  - (Iqf,2  + X1X)'1  w2  . 

It  is  well  known  that  an  X satisfying  (7)  exists  and  is  unique  if 
and  only  if  A^  and  A0  have  no  eigenvalues  in  common.  In  practice  only 
such  cases  interest  us  but  we  want  the  algorithm  to  be  robust  in  the  face 
of  some  perverse  or  extreme  requests.  Clearly  if  A1  = A2  we  want  the 
algorithm  to  do  nothing  rather  than  to  fail.  In  such  a case  s = 0 
which  is  far  from  Invertible.  By  taking  £ * 0 and  setting  C^/f,  * C0/.f.  = 
* X = I the  algorithm  will  work.  When  the  eigenvalues  of  A1  and  Ap 
are  close,  in  some  sense,  then  £ will  be  chosen  so  that 

max{£,IXI)  - 1 approximately. 

There  are  infinitely  many  C's  satisfying  (8)  and  (9)  and  any  of  them 
will  do.  In  the  absence  of  other  constraints  the  symmetric  solutions  are 
the  natural  ones;  if  cj  * , C2  =■  C2  then  = Sp,  but  this  fact  is 

not  obvious.  In  this  algorithm,  however,  we  prefer  to  choose  Cj  and  C„ 
so  that  A1  and  A2  have  a convenient  form  for  most  applications. 
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It  is  not  necessary  to  compute  and  Sp  explicitly.  Write 


e,  - 


’-Sl  C2 

s 

' l2  0 ' 

1 

X 

C1  S2 

o ^1 

l ci  x J 

^ - C^/C,  the  scaled  version  of  C^.  Then 
(10)  P * 


and  P is  best  applied  in  this  factored  form.  In  practice  the  orthogonality 
of  P is  completely  determined  by  the  accuracy  with  which  the  C's  satisfy 
(8)  and  (9). 

It  is  not  necessary  to  compute  A^,  A.,,  or  B explicitly  since  they 
will  emerge  when  the  similarity  transformation 


(ID 


”1 
0 A 


is  effected.  For  completeness  we  give  the  formulas 


s,  - xl\ 


S2  = 


(12)  A,  - ^A^j1  , A l = CjAJC"1 


B = A2C2TsJ-SlCl  A1 


'"l-*  _ c^a  r"l 

? "2^2  111 


T 


r 


8 

[A,  B' 

4.  The  Algorithm  for  SWAP  A , 

~ [ 0 A2 

1.  Clear  the  (2,1)  block  of  A. 

2.  Solve  AjX-XA^  ■ f.B  for  X and  £ using  subroutine  TXMXT. 
i.  is  chosen  so  that  I XI  4 1 

3.  If  C ■ 0 then  exit. 

4.  Solve  Ejc,  - (cWf1  = W1  for  using  CTCEOW. 

5.  Solve  C]l2  - UWX)-'  = W2  for  C?  usinq  CTCEQW. 

6.  Premultiply  A by  P usinq  NEWCOL. 

7.  Postmultiply  PA  by  PT  usinq  NEWROW. 

8.  Update  the  matrix  of  orthogonal  transformations  using  NtWROW. 

A A 

9.  Force  the  diagonal  elements  in  the  new  blocks  A1  and  A0  to 
be  equal. 


Name 

Executable  Statements 

Count  for  2 * 2 Case 

SWAP 

17 

32n  multiplications 

TXMXT 

52 

42  multiplications 

CTCEQW 

17 

32  multiplications, 
4 square  roots 

NEWCOL 

22 

16  multiplications 
per  column 

NEWROW 

22 

16  multiplications 
per  row 

i 


When  A 


‘1  i 
'l  ^ 


^ 5.  Solving  AjX  - XA,,  * B 

1 * 1,2,  the  linear  equations  which  determine  X 


can  be  solved  stably  in  closed  form.  Let  6 » c^-o^,  then  the  equations 
may  be  written  as 


9 


i 


(D 


Yl!2 


_ 

X11 

bll 

Bl!2 

x = b ; x = 

x12 

, b = 

b12 

C 

X21 

C\J 

jO 

X22 

l b22 

where 


(2) 


'6  -Y2 

r2 

*^+B2Y2 

-26y2 

C = 

. -b2  6 - 

, c = 

-26B2 

62+62y2 

Multiply  (1)  as  indicated  in  order  to  make  the  coefficient  matrix  block 
diagonal , 

, ,2 


(3) 


c‘-Vi  0 L _ f c -Bi  1 


0 


x = 


*i  c 


b . 


Now  let 


6 = (C2-B]y1 )_1 


t 26>2  ) 


2662  t 


/d 


where 


(4) 


T = 62  + 82y2  - 6-|Y-,  , d = T2  - (2<SB2) (26y2) 


and  premultiply  (3)  diag(G,G)  to  find 


' G 

0 ' 

C 

V 

. o 

G 

■ 'Y1 

C 

(5) 


’ g :o 
. o ;g 


6b 
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Y2b12  " Plb21 


"B2bl 1 + 


6b 


12 


" 6lb22 


"Ylbl 1 


" Ylb12 


+ 

- Bob 


6b2i  - Y2b22 


2U21  + lSb22 


♦ (26  -T)>2bi2  “ 

(262-x)82bjj  + 4>6bj2' 

-T>1bn-  26ylY2b12  + 


T0lb21 

26P]82b?1 


26> 28,b2? 


Tpl b22 


(J>6b91  + (26£-T)Y2b22 


TYlb12  + ^2(S  -1 ) p2b21  + 


26Y162b1 


defining  y,  where 


4>  = t - 2>202  = 5 " (Pi>i  + ^2>2^  % 0 ’ 

<P  * 262  - t = 62  ♦ (61y1-62y2)  • 

Inevitably  (5)  is  Cramer's  rule  and  d = det(A^®I  -I®A0)  so  that  d = 0 

if  and  only  if  = a2»  8^Yi  = 62Y2- 

Among  all  the  coefficients  in  the  linear  combinations  of  the  elements 

of  B which  are  given  above  only  t and  involve  genuine  subtractions 

and  possible  loss  of  information  through  cancellation.  However  by  rewriting 

them  in  a more  complicated  form  all  unnecessary  loss  can  be  avoided.  From 
2 2 

(4)  x = 6 +P2Y2”81Y1  an<*  ^ eitber  of  t<!  or  -B-j >i  is  tiny  compared 
with  the  other  two  terms  we  want  to  add  it  in  last.  Similarly  for  i|i.  Thus 
we  use 

2 2 

<J>  = (81y1  + max{6  ,-B2y2))  + min{6  ,-62y21  , 

(6)  2 ? 

t = (82Y2  + max(6  ,-B^))  + min{6  ,-t^^ > . 

Here  is  an  example  for  a machine  with  a relative  precision  of  8 decimals, 
i.e.  the  floating  point  result  f 1 ( 1 0®-9)  is  108  whereas  fl(108-10)  is 
1 08  - 1 0 = 10(107-1).  Let  62  = 9,  6^  ■ -(108-10),  82Y2  - -108  then 

from  (4),  computing  from  the  left,  t = fl (fl (9-1 08 ) + 108  - 10)  = -10  , 

from  (6),  computing  from  the  left,  t * fl (fl (-108 + (108-10)  + 9)  = -1  . 
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I 


KT 


warn 


If  we  are  given  a matrix  M with  eigenvalues  near  «1031  and  are  evaluatinn 
exp(lOM)  then  values  like  the  ones  given  above  will  occur. 


Norma  Hz ation 

The  important  matrix  in  effecting  the  orthogonal  transformations  is 


' l2  0 

UJ 

h- 

X 

1 

■ 0 

( c x J 

and  we  want  our  formulas  to  be  accurate  right  out  to  both  extremes: 


’-10 

0 1 

and 

0 I 

I 0 

An  appropriate  way  to  achieve  this  is  to  choose  E so  that 


max{£,IXI } 4 1 . 


Equation  (6)  above  yields  y so  that  the  correspondino  2x2  matrix  Y 
satisfies 


Aj Y - YA?  - dB 


where  d 

is  given  in 

(4). 

To  get  x 

and  E 

Case 

J.:  n £ d. 

take 

x ■ y/d. 

C ■ 1. 

Case 

2\  n > d. 

take 

x * y/n, 

E = d/n 

«yl„,  then 


The  Algorithm  for  TXMXT 

We  solve  A^X-XA^  ■ £B  when  A^  and  A?  are  2*2  standardized 
matrices  as  follows: 

2 

6 ■ cij  - • 6sq  ■ & , 

"1  * Vl  • n2  B2y2  • 

e ■ <f><S  * <S ( 6sq  - (ttj+Tt^)) 
f1  * e (n1  + maxlfisq.-Ttgl)  + minUSsq.-i^l  , 
f2  ■ t * (tr2  +max{6sq,-n1 })  + minftfsq,-^ ) , 
q - 2Sy2  , h - 2<5B2  . 
d ■ f2  - qh  . 

At  this  point  y can  be  evaluated  from  (5).  Then 


n - lyl^  . 
x - Cy/max(d,n)  , 
new  £ ■ £*d/max(d,n)  . 
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6.  Sol  vino  C C = W 


In  some  applications  the  A^,  1 « 1,?,  have  the  special 


form 


N ■ V?  * 


0 B, 


3 .-y . < 0 , 
11 


and  we  want  to  have  the  same  standardized  form  (equal  diagonal  elements) 
Because  A^  is  similar  to  A^  we  must  have 


fli  ’ “l'2  + 


foe,) 
9,  o 


Vi  * 'Vi 


This  requirement  fixes  the  matrices  and  of  the  previous  section. 
A straightforward  way  to  derive  formulas  for  and  is  to  obtain  a 
particular  solution  to  (8)  via  the  Choleski  decomposition  and  then  to 
standardize  the  resulting  diagonal  blocks. 

Let  R.|  and  R^  be  upper  triangular  and  satisfy 

r|r.|  = = U2i2  + xxt)_1 

R2R2  = W2  = (6ZI2*XTX)"1  , 


where  X solves  (7),  A^X-XA2  = f.B.  Next  defi 


ne 


a,  i 


A2  ~ R2  ApRr)  . 


Now  let  and  J2  be  the  unique  plane  rotation  matrices  which  stan- 
dardize A^  and  A^t  l.e.  both 

*1  = and  A?  = J2A2J2 

have  equal  diagonal  elements.  The  appropriate  C1  and  C?  are  therefore 

^ * cyr.  = - c 2/r.  = J2R2- 


Let  us  drop  the  subscript  and  dot  from  and  A y The  condition 


CTC  =■  W = (f,2I2  + XXT)_1 

imposes  three  quadratic  relations  on  the  four  elements  of  C.  If 

’ a B 1 

A * y a then  requirement  that  CAC  have  equal  diagonal 
elements  (both  a)  Imposes  another  quadratic  constraint,  namely 

Pcllc21  * yc12c22  ’ 

which  suffices  to  determine  C.  However  the  direct  solution  of  these 
nonlinear  equations  is  far  from  obvious.  Instead  we  shall  derive  the 
solution  in  a straightforward  but  lengthy  manner  via  the  Choleski  factori- 
zation of  W.  The  final  algorithm  is  however  very  compact.  Let 

d2  = detU2I2+XXT)  = S4  + C2(HlxiJ|2)  + (det  X)2  . 

Then  define  M by 

u m/h2  - l[c2  + X21+X22  '(xllx21+x12x22) 

W * M/d  = -W  222 

d _ “(xi i x2l+xl 2X22 ^ + X1 1 + x12 

and  note  that 

R b 1 ’'Si  m12/’STl 

d[  0 d/vfiyj' 

2 

is  the  Choleski  factor  of  W.  Note  that  det  M 3 d . The  next  step  is  to 


form 


ri!  + aI2  • 


Now  let  J be  the  plane  rotation  which  standardizes  A 


6(c^-s^)  - 2Ssc  2fisc  + a.|2c?  - a?1 

26sc  - a'1?s2  + a21c2  -<S(c2-s2)  + 28sc 

where 


a - (a12+a21)/2  . 


The  proper  choice  of  c * cos  G is  therefore  given  by 


So 


tan  26  = 2sc/(c?-s2)  = <S/a  . 

2c  ^ = 1 + cos  2o  = 1 + a/\>  , 
v * A3  + aT  , 

c = /0  +Tal/\'T/?  , to  keep  |0|  < n/ 
s 11  sin  20/2  cos  0 ■ <S  sign(S)/2cv  . 


Our  object  now  Is  to  get  rid  of  the  Intermediate  quantities  and 
express  C in  terms  of  d and  M.  So 


« * (0rH  ’ Y(rl2'r22I)/2rllr22  ’ 

- yld*"  - (m2?  - t^1/y)|/2m1  ^i  , 

= YC/m^  , defininq  r., 

<S  * Vr^/r^  “ » 

v * IvU/nt^  where  <t>  ^ /f7~+  m^0  . 

Since  fly  < 0 the  expression 

2 2 2 
u)  = m.p  - itm. . /y 

is  positive. 

At  the  cost  of  an  extra  square  root  the  important  quantity  <\  can  he 
written  in  a form  which  is  attractive  for  finite  precision  computation 

C s [d2  - (m‘2  - Hm^/Y)l/2d  - (d-u>)(dM/?d  . 

Having  computed  d,  M,  C,  and  $ we  obtain  the  desired  formulas: 

o * m^/2d2  , 

c^  ■ cr^  ■ f C | /<£}  , 

p 

C21  ’ srll  " rll  sign(a)/2v(cr11 ) * siqn(f  1 . 

c12  " cr12  " Sr22  " ^cllm12'c21d^/m11  * 

■ <c2il"i2*cnd)/",n  ■ 


c22  ' sr12*cr22 


V 


For  completeness  we  note  that 


0 

ft  ' 

a 

C11  C12 

OB’ 

C22  ~c12 

Y 

0 

C21  c22 

. Y 0 

• ~c21  C11 

so  that 


ft  3 (flc^-Yc^2)d 

Y * fiY/ft  • 


The  matrix  C is  computed  by  the  subprogram  named  CTCEQW  (i.e.  C C = W) 


Computation  of  C 2 

The  subprogram  which  computes  Cj  from  d,  X,  B,  y can  also  be  used 
to  compute  C^.  Recall  from  (9)  that 

cjc2  - (1 2 + XTX)_1  . 

By  symmetry  d?  * det(I  + XTX)  = det(I  + XXT).  Moreover,  from  (12) 

A2  * . 

By  transposing  the  data  we  can  use  the  same  formulas  as  qiven  above  for 
The  data  is  d,  x\  Y2»  ^ and  the  outPut  be  C2,  y2»  ln  other 
words  it  is  only  the  interpretation  of  the  parameters  which  distinguishes 
the  computation  of  C2  from  that  of 


7.  Performing  the  Similarity  Transformations 
In  practice  and  A2  will  be  contiguous  submatrices  on  the 
diagonal  of  some  big  block  triangular  matrix.  The  similarity  transformation 
determined  by  P affects  elements  in  the  same  row  or  column  as  those  of 
A1  and  A2  as  indicated  in  the  figure. 


£2(Cv-XTu) 

^(Cu  + Xv) 

Notice  that  the  number  of  multiplications  required  to  effect  this  is  pq 
for  each  of  Xu  and  Xv  plus  q and  p for  the  application  of 
and  C.j.  This  is  the  same  as  for  multiplication  by  the  full,  non-factored 
version  of  P except  for  the  (p+q)  multiplications  involving  f.. 

There  is  a surprising  difficulty  in  writinq  a program  to  effect  this. 
The  program  must  work  for  any  values  of  p and  q and  this  condition 
prevents  us  from  supplying  the  input  data  as  values;  they  must  be  names  or 
references  since  the  number  of  them,  p+q,  is  not  known  at  compile  time. 
In  other  words  the  subprogram  is  informed  that  elements  m+1  through 
m+p+q  of  an  array  Y are  to  be  transformed. 

The  disadvantage  of  this  constraint  is  that  the  same  code  cannot  be 

T 

used  for  effecting  the  postmultiplication  by  P . More  precisely,  the 
price  of  using  the  same  code  for  both  cases  is  a loss  in  elegance  and 
efficiency.  The  difficulty  can  be  seen  clearly  by  looking  at  the  listinqs 
of  the  subprograms  NEWCOL  and  NEWROW.  They  differ  only  where  a variable 
Y(i,k)  in  NEWCOL  corresponds  to  a variable  Y[k,il  in  NEWROW. 
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8.  Gaussian  Elimination  for  Solving  A^X-XA^  - B 
The  linear  equations  defining  X can  also  be  solved  by  block  Gaussian 
elimination  in  about  half  the  time  required  by  the  algorithm  just  described. 
Three  different  factorizations  are  appropriate  (i.e.  stable). 


Case  1:  6Z  » max(-0^Y-|  »-62y2) 


I2  0 C 0^2  ^ 

Y-j C I2  , 0 C-B^Y^C  k Xj 


Case  2:  |y  I > 1^1  » max(6Z,-B2Y2) 


’ bl 


r2  0 

Cy”1  y’1 
^Y1  Y1 


0 -(C-BlYl)  Jl 


Case  3:  | y2 I 1 I B2 I » max(6Z,-B1Y1 ) 


-1  -1 


0 C2-82y2 


Cl  x 


xn  j = f X12  1 g = bn  g 

X12  ’ ~2  x22  ’ b21  * °2 


In  each  case  X can  be  found  with  16  multiplications  and  4 divisions.  Further 
rearrangements  should  be  made  when  | S-j  | > Iy-j  I in  Case  2,  |82l  > |y2I 

In  Case  3. 

The  extra  length  of  the  code  (100  statements  versus  50)  does  not 
appear  to  warrant  a saving  of  16  multiplications. 
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9.  Swapping  Large  Blocks 

The  algorithm  we  developed  for  swapping  was  quite  general  with  A1  p*p 
and  qxq.  However  the  individual  subroutines  TXMXT,  CTCEQW,  NEWCOL,  and 
NEWROW  were  specialized  for  p < 2,  q < 2.  Here  we  want  to  point  out  that 
general  versions  of  these  programs  are  readily  produced. 

1.  A^X-XA^  = B can  be  solved  for  X by  the  algorithm  of  Bartels 
and  Stewart  [B  and  S,  1971].  In  our  case  A^  and  A^  are  already  in  real 
Schur  form  and  X can  be  partitioned  to  match  A^  and  A^.  If  the  equa- 
tions defining  X are  taken  in  the  proper  order  the  system  is  trianqular 
and  can  be  solved  by 


A^X.,  -X.  = 

kfc  kfc  kt  u 


kl 


f a^x  + \ x a!^ 

j=k+!  kJ  J*  i=l  kl 


The  proper  order  is  k = p,p-l,...,l;  l = l,2,...,q.  Here  p and  q are 
the  block  orders  of  A1  and  Ag. 

2.  CTC  = (£^+XXT)”\  The  positive  definite  matrix  t}  + XXT  can  be 
formed  explicitly  and  its  Choleski  factorization  R^R  computed  in  a standard 
manner.  Then  RT  can  be  overwritten  with  its  inverse  to  qive  a solution  t. 

3.  The  execution  of  the  orthogonal  similarity  transformation,  in 
factored  form 


no 

o 

’ -XT  c ' 

l o e2  ] 

, 5 X 

presents  no  difficulties. 

We  mention  this  possibility  only  to  reject  It.  The  rival  method  is 
simply  to  swap  A1  and  A^  subblock  by  subblock,  using  the  proqrams  which 
we  have  presented  here,  that  is  by  swapping  many  1 *l's  and  2*2's.  The 


p 

operation  count  for  each  method  is  approximately  (p+q)  n multiplications 
and  additions  but  the  general  procedure  sketched  above  would  require  signi- 
ficantly more  program  statements. 

In  the  language  of  computer  science  we  are  recommending  the  recursive 
swapping  of  big  blocks. 
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1 0 . Test  Results 

(a)  6*6  Matrix  (Separated  Eigenvalues)* 

1.  Original  Matrix 


2.0000 

3.0000 

4.0000 

5.0000 

6.0000 

7.0000“ 

1 .0000 

2.0000 

5.0000 

6.0000 

7.0000 

8.0000 

6.0000 

7.0000 

8.0000 

9.0000 

8.0000 

9.0000 

10.0000 

12.0000 

11.0000 

-1 .0000 

12.0000 

2.  Swap  1st  and  2nd  blocks,  2*1  case 


"6.0000  -4.2583  -4.6036  9.6667 

11.328 

12.990  “ 

2 . 2 - 1 0" 1 4 2.0000  3.2930  -3.8212 

-4.3070 

-4.7929 

-0.91103  2.0000  -1.3978 

-1 .4569 

-1 .5161 

8.0000 

9.0000 

10.0000 

12.0000 

11.0000 

- 

-1.0000 

12.0000 

3.  Swap  3rd  and  4th  blocks,  1 *2  case 

f 6.0000  -4.2583  -4.6036  13.192 

-13.265 

6.3656  H 

2.2 

* 10"14  2.0000  3.2930  -4.8713 

5.1785 

-2.3615  ' 

-0.91103  2.0000  -1.5437 

1.8492 

-0.75551 

12.0000 

0.69449 

5.5837 

-15.839 

12.000 

-4.5244  ' 

4.2  * 10‘14 

-1.3*  10' 

14  8.0000 

Computations  oerformed  on  14  digit  machine,  results  rounded  to  5 figures 
for  display. 


4.  Swap  2nd  and  3rd  blocks,  2*2  case 


6.0000 
2.2  x 10' 


13.402 

-14.062 

-1.3625 

-3.1781 

6.3656 

'14  12.000 

0.63866 

-3.6006 

0.39991 

5.3584 

-17.224 

12.000 

-0.21201 

0.40812 

-5.1223 

-1.0*  10"13 

4.6  x 10'14 

2.0000 

2.7985 

-1 .6498 

1 

-4.7  x 10'14 

2.5x  10"14 

-1.0720 

2.0000 

-0.35352  ! 

4.2  x 10" 

14  , , 

-1  .3  x 

rr 

i 

o 

8.0000  ! 

— J 

(b)  6x6  Matrix  (Close 

Eigenvalues)* 

nal  Matrix 

6.0000  10"4 

4.0000  5.0000 

6.0000 

7.0000  n 

1.0000  6.0000 

5.0000  6.0000 

7.0000 

8.0000 

6.0000  7.0000 

8.0000 

9.0000 

6.0001 

9.0000 

10.000 

6.0001 

10"4 

-1.0000 

6.0001 

2.  Swap  1st  and  2nd  blocks,  2*1  case 


6.0000  0.99984 

-4.9995 

-5.9992 

-6.9990 

-7.9989 

6.0000 

4.0006 

5.0010 

6.0011 

7.0013 

-2.4996 x 10'5 

6.0000 

7.0000 

8.0000 

9.0000 

6.0001 

9.0000 

10.000 

6.0001 

10"4 

-1.0000 

6.0001 

♦ 

Computations  performed  on  14  diqit  machine,  results  rounded  to  5 fiqures 
for  display. 
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3.  Swap  3rd  and  4th  blocks,  1 *2  case 


6.0000 

0.99984 

-4.9995  - 7.9996 

5.9992 

6.9983 

6.0000 

4.0006  7.0019 

-5.0010 

-6.0004 

-2.4996 x 10'5 

6.0000  9.0008 

-7.0000 

-7.9991 

6.0001 

9.9992  x 10'6 

0.9999? 

-10.001 

6.0001 

8.9991 

- 

CC 

1 

o 

X 

-4.3  x 10"19 

6.0001 

4.  Swap  2nd  and  3rd  blocks,  2x2  case 

"6.0000 

-4.9997 

-0.99972  -5.9991 

-7.9995 

6.9983 

6.0001 

2.4995  x 1 0'5  7.0000 

9.0008 

-7.9992 

-4.0008 

6.0001  -5.0011 

-7.0020 

6.0006 

-2.1  x 10'23 

-6.6  xlO'24  6.0000 

10.001 

-8.9989 

8.9 x 10'24 

-3.8  xlO"25  -9.9994  x 

10"6  6.0000 

1 . 0000 

3.9  x 10‘ 

18  -4.3  x 10 ~^9 

6.0001 
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11.  Program  List  Inc 


*UHGiUi“l  M'  C*AI>(NW,N,7,  ’ , J1  ,11  |L<’  ) 
n I .1  : f.  - ( \'*,  N 1 , " (*n* , M » 

KE.*L 

rOU!«M.  Nf!  < xr  | , t > » v<  1 « l>  > •(  C 1 ( l « I > , v<  I ( 1 > ) « < C,  < t . I > « v ( 1 , . ) I 


PXrHANf,1  UJJAl  ‘ N T kIACONAL  f'lffK'  Tl  A 
rriTw  JI'«  ltl'l’Y  T ’AM,,nMVAMrNrl 

Wf.m  v*.  r '"jm  i'»  ▼ . "i  .*k  ti  T . i l 

l « j 1 ♦ i 1 

t?  .n  - » 

z "■  l.J 


(•  (,  INN  INC 
.MO  IN  l>  , 

, i ir 


If  W**  J1  OY 
nr:;  v I N i TH 

-IV  t.  ? • 


r v*»*****xi*****(.'Lr/»K  thp 

L„  OC.  .5-JxJUL  1 

D P Mltl? 

T ( I > ♦ T,  •!  »J-  I IrC, 

0 

r **  * * ■*  A I CK 


m oc> . 


CAI  l.  TXMXT  CN*tN,T,Jl,  J3,l  l,t 
!*•  ( * ,r‘0.  % 1 pr  TlJfM 


>_X  . i;  « T 1 , XhC->r 


7 , X ) 


^-lAti^A^dAicrt'Mn.rr  n XHTOf  Cl'' ’'Cl  * (70 O'*  l * X * X T 1 

C A*  i'  \v .« v>  > a m n ^ * • ( i ' C .1  y v r ; » ( f o 1 I ♦ x r * n \ * - 1 , 

C All  -’Cr  ow  » ',  x ( i , i 1 , x ( I » . x(  I , ? 1 , v(  a,  p)  , ’(  1 1 , j?»  , 7(  l ’ , J | > ,C  I ,t  I ) 

r hi  ' rr-  , i ) . x ( l , ' ) ,x  < ? , 1 1 , v ( * , ? 1 , T f ,|4  , j ? l , T ( i r,  ) a ) , c i ■> ) 

r 

' C*****V»  **  **  * «i*pp  »?»-OKM  T o ANSrorfMAT  ION  ON  COLUMN*?  ANf>  »0*1  I"  7, 

C * * * ■*  * * ■niAiA  ti  + yan  a ff  P. 

f Al  I N'xvrr  (7,1  1,L?»Y,T,N",N,  11,  1 ) 

CiLLL  vr  IYLC  UlU  ,tr,y  » J, »A  , Ji  f NW> 

CAlt  N ' w V f C (7,1  l, l >,Y  ,*  ,NM  ,N  , Jt  ,NM) 

c 

C •***<  »'  iv.T  «*,  p yi  I otjA  l T r Y ’(  oi  AC.rNAI  • LI  Mf  NT  % IN  FU  OCK.  ' . 

ir  (LJ  iio,  i»  ) r ( J 1 . J l >-  ~(  J 1 ♦ 1 , J 1 ♦ 1 » *(  T<  J 1 » .It  I ♦ T(  J !♦  1 , J 1 ♦ l 1 1 /<', 

T r Tn  .r  o.  a ) t ( ja  -i  . ia  -i  > *y  < , j*  >r  (t  < ja-  1 , J4-|  ) ♦ (.1  A,  J A ) ) /?. 


OPTION 

FM* 


'•UMWOUTIM*  ( >,X??,l1F7A,r»M,C,L) 

DIM1-  f 4 C I C N c (.»♦?> 

c'^TniT  aKi  An6f»n{TSlTT«!  SAlUY TTV  f f'rt  CY*r  • x • (/so*i  ♦ K*>xt>**-i 
71  c«  7*  7 

ir  (l  .r, T.  |i  c.n  rr  to 

C<  ^ )«  i ./con  y ( 7?p ♦ x u * x 1 1 ♦ * ? \ ! x .*  1 ♦ x i ?* x i ? » 

1 J r-  M l I « 7 f.n  ♦ X ’ l * * .7  ♦ x .*  • r> 

• vi  .>  * - (x  1 l *■  x l »x  l x > 

0«LM1  y i i v * >♦'1^  •*;>  ) tI:‘«W.**2 

'G  A ■ (fMi;1**.1  - or  r A»v  M 1 1 +*  ?/  f.AM  I 

*r  ta  • (pro  - rcAIMPTO  ♦ r GAI  /<?,  -nitiTiM 

"fc  1« 5 JL'ULIf 

► a r »r  m i | / f * . fw  n) 

>*(  1 . 1 ) = SOf-  T ( I AC  M 1 . .'lA'ISl  7 7 » | /PHI  » * 

C ( s’  . I l-c!r>*'||.',?rTMM'M.'*IAI,r|l'MIV(  | , | ) ) 

C ( 1 » - ) ;Ut  1.J.I  )*_NiZ-CL^.n  ■ T,')7’’V|  i 

c(f . .7 1*  ( r ( >,  i »**-m  | .*♦(-<  i , n *k  rni  /f  m 1 1 

W T1)0N 

r no 
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MU-PTIT  I NT  TXMX’  (My  ,N,  T , X) 

O r V,  N\,  triN  T ( NM,  N)  , P ( NM,  N ) , * ( 2,  ? ) 

C SOLVE  P0«  LI  MV  t.?  MATRIX  x IN  T1»X  - x*T?  = /w, 
r t i \no  t?  o«-r, in  in  ji  ano  j? , 7 c;rv‘N  on  *ntpv 

C 1UT  °\  |-  X !T  r IS  CHA^G^rt  rn  FNC(JPF  noc  m ( y ) J t « , J , 

miai  =x  t > .1  ) -x  (i ,?  ) = \.(J242)so . 

IF  ( T . -*0.  f . > wtr  T » jr^ 
nr|i'(  1 1 . J 1 ) -T  ( J ?,  .1  •> ) 

Kr?*Ll»l  P -'■* 

LI — LQ-J a jlXt  - 1)  jL  TC  r 

r 

C*******«******TJ  ANO  Tp  M4  Vr  TH<*  SAME  f I Gf  N VALl/'-'j  , KCTUPN  /=3. 
« X I l * l > = X<  3,?Jr  l. 

c 

c «'*v  *>:-•<*■*  x>  r.  MI  Nr  DIMI\STCN^  DF  f CLUT  I TN  X. 

_=> 50.  T J -L-LC-».2C.«  5C  «.ACX».K 

C 

f*.K.F*«***«***r  l IS  1 PV  t«  T 2 I'  1 riv  1. 

JO  x MAX = APS  (R(J1,J?)I 

n — ar  st  nr  i \ 

x ( i , 1 ):P(  J),  j?)*(  x/am  a xt  ( x -4  a y,  )»  )<  p i r,N(  t o ,n:  i I 
.()  T 1 TO 


O ‘ A-*  ''  X V F.’X  ''■F*T  I T S I MY  1 


T2  IS  2 ilY 
3 0 0*0n_**?-T<  J?  , J2M  )*T(  J2  A I ,J  >) 

X<  1,  l >*<p<  JI  , J?)*0FL>O(  Jl  , J3H  ) * T(  J?f  j , jjl) 
XU  * < P ( Jl  . ) *T  ( J?  , J?  ♦ J.  > »BI  ( J 1 ^ J 3 M > *Pf.|  ) 

GL  th  sc 


C lfr  trt  *c  in  s«  A *.  + « * • t a *v  y 
UJ-C-l 


I IS  ? M Y 2 1 T.J  T < 1 r,  Y l, 

( Ji.i  UJi  ) 

x(  1,1)  *(  0rL<tr  ( J l , J2)  -T(  JI  ,J1MI*«UIM  , JP  ) ) 

x ( i j - ( -t  r iim*ji>*‘3<ji»j?>  ♦hel  wuin.irn 

XMAXsAMAX!  ( A|lc  ( X f 1 ,1  > > , APS (X  < 3 , 1 )>) 

-GL._ia.sji  — ..  . 


a r 


ii  i’  x # «***«  ▼ j |c,  ^ uv  2,  T 2 I i •>  ,<Y  ?, 

PI.  T 1 *T  ( j i , j j ♦ j ) 

mtY-1  xIIJl  »UJD-- 

PF  T 2“  T(J?»J2A1  ) 

GAVP*T<  J341  , J3  ) 

P 1 *PET1*GAM  1 

B2 F.B.ET  P ft  f. AM2 . . 

r>  <^n=‘')r  i * *•  p 
l-.r'LM  0~<J-  (n|*n?)| 

ri  = ( PI  ♦ AM  AX  1 (Dp-I),.n?  ) ) ♦ AM  IN  1 ( OSQ  , -P  ? > 

£-2  ?XP 2 tAVAXIlUiO  J_1  (Q50  , -PI  ) 

H-2.0  *0rL  »'irT2 
r,  r 3.  0*P.;  I >■  r,  AMP 
,A  -t  >>:•>:  0-  .-,■*>  M 

! I ( P *F  I • C • ) CC  ’■p  A 
‘ l 1 *P(  JI  , JP  ) 

- 1 P-P ( JI  , J ♦ I > 

«21  *P(  J 1 M t J2~> 

P2?aP(J14 I , J2+ 1 > 

x ( i « i lip  i i *p  ♦«  t ?»ga  i-rxiTn^Tixr  ?-pj,’*Br*i*G 

X ( l , ? ) =Pl  l A nrT?Tfi  4PI  ? ** .pi  I * Pf-T  1 *H-P  PP*'  PFT  1 *r  ? 


x ( :T."fV  -f  : I ''  >*,/.  m i>'  i i ; * c» a 4 j"* rt a f.  j-f  ♦ f 
x ( p , p ) - j i *<.;am  | |?»G\>i|pr?»F  ■>  j < r,c 

XM,\  X ; A MAX  1 ( AfP.  ( X ( 1 ,1  ) ) , ( X ( J , > ),/'>'• 

r 

OR'^  * *'  «WW9++ * A « * r xi  si  j fir  nopmix  ) < l* 


< 


<•• « , j 

I ♦- 

’ . t ) ) , A 


( * ( 


*4  P 


f •) 
T ■' 


fee  7 /AMA  X | ( XMA  X ,P  ) 

on  «•<>  j*  i , t 2 

j>_a  *o  i f i » y i . _ . 

x < » , n - x ( » , j » 

IF  (X^AXoVoD)  »:/^/«MJx 
P1-  TUPM 

.ifclC 


*-'.UPc,r')'i  wr  Nr^vrc  if  ,m  ’,v,nw,n,vi  , iMf  i 
OlMf-M.ION  Y{  I * , 7(  4)  , M A > 

r 

c c qmp U t*  :»■'  i r*  v - = MifcU.  •ci.M£*-u..  t..  a^vi  * nlwv 

C KHMr  0 1 IN  N 1 BY  Ml,r?  IS  N 2 n Y N ?, ANO  I U,  V)  *(  Y(  M I ) , Y(  Ml  4 l 1 , o » « > . 
c THr  MATrtCf  s Y,Cl,r?  A^T  STOPfn  tN  7,  F IS  SCAI.  INC  F AC  Tf  o Fit'’  X* 
r 4‘«''  N 'N'tl, 

c m*'n»  ml oc*  ‘-xOHA'ir.*'  T*'  ANSFr,FM*T  IQ1'1  rN  oiiomns  m | to  n of-  v 

• IN  THF  'NlFNi't-nw'  <‘TArTI\r,  »PH  mi. 

C NY  = t'*l-l,Ml*  •-  '*  1 - 1 ♦<  M 1 - 1 1 *N'«  = (V  1 - 1 I -t  ( N«-*4  1 | , | l':V)  , .1  = N M, 

r *H  N IKFjNM, 

C PLr-  f-'QRM  i.tL'X<  L-XC*i4\GL  T‘  A N SI" U fc". V A. T LC N LN  I 1 • - r N f -W'J  . rJI  > 

C IN  THF.  Ml  4N7  Cijt  L'M4  S NT  ANT  INC,  » IT  h OH'JMN  Ml, 

C NY  * ( = 1 ♦ I Ml  -■»  NM  = ( M|  -1  1 * N*41  -KM,  l h»l  , J*l« 

1=  TNT 

_ lq»  j-nm-t  ♦ i 

II  (I  , l * . J 1 l P e M 1 
NY  r (MJ  - 1 1 M»*  «l  i - » 

l stv  1 ♦ ’<  MO- 2 

GC  T0  ( ld«  2C  ».  JC  • A LI*  l — 

c 

C***tiN**^Mvr»Nl»l|  N 2 * l • 

1)  or  IS  *C  * L O « N 

U L}JFJtlNYl2»Il*PsZ.U  »1  >*Y  (NY  tl  ) 

»(?)iY(NYH»*F  4 J I 1,1  * V Y ( N Y 4 ?'•  I 1 

V ( NY  f 11  s/(  I ,M"WI!  I 
VIM  v*2'  r 1 = 71  I . *>*» .12) 

-iJ liX-J^CJU. ......  - 

wr  TUMN 

c. 

C**************Nt  *2  , N2»t« 

20  DO  23  K»LB  .N - _ 

W I \ ) = Y (NYT.T*n*r-7(l,l  1 *Y(  NY  ♦ I 1-7  ( 2,  1)*Y(MYT?*1  > 
tel  ’I  =Y(NYTI)*FTZ(  1 ,1  )*Y  (NY  AIM  ) 

W I U = Y(  M Y ♦?*  I ),,F-47(2*l)*Y(  NY»  3*  I * 

t-LI^LUlr  ill  .:  ) v,»  ( LI 

YINY*?«1  1 = 7 ( 1,3  >**(2142(1 

V (NYM'-I  1 = 7 ( /,  .3  1*4  ( 2 > ♦ 7 ( 2 , A ) *w  ( 3 > 

?s  N Y =N  Y ♦ J 

rLL-UON  . ..  ...  . 

r 

C«»t  ,*44Y»*  l cl  , N 2 = ^ • 

30  nr  3*s 

HI  1)  »Y(NY42«I  )+r_-U  1 .1  ) T y ( NY  4 j j 

W(2>«V(NY4  3Ml*l'-/(  l,?t"V(NY4l) 

*(3  » »Y  (NY ♦ I»*r ♦?(  I , l )4Y(NY4?»t  )4’(  1,  5|'Y(NY4I*I1 
YIN  Y4-  I ) « 7 ( 1 . r 1 * *(  l 1 4 2 ( 1 , I*  * ( 2 > 

V ( NY  tZ"  1 ) 1 Zl?^  ) »>g.l  > »ZJ  2.  £ ) **1  2 \ 

YI'jY4  .«i|  1=7(1  , T)  **(  3 » 

15  t Y:'IY4  ) 

r.  r -y  . il  n 

C*  ♦"»  * * * **  Y 5 «iMT*N  1 • r , N?«2* 

4 C or  4 c K*l.P,N 

*1  1 ) = Y(  N Y4  3*.  1 ) *►--  7 ( 1 , 1 1 * Y ( NY  4 I ) - 7(  2 , 1 1 * V (MY  42*  1 ) 

I*  ( 2 1 = Y(NY  44  * l ) *r  . 71  1 , >)*Y(NY4|  >-2(  ?,  7)  " YIN  Y4?»  1 ) 

M31  =Y  I NY  ♦ 1 1*4"  4-2(1  , 1 1*Y(\Y  4 V>  I >4  7(1  ,7)«Y(NY444l  ) 

. _ vtJLf-  ) =_Y L'r  V+Z*  urr  *U  i • 1 > ,vYCNYtJ*  ! > ♦ < J,..'  )*  y (nt***  I > 

V( MY  411=7(1, * >*W(1  >4/(  1,4  )*W(2> 

YIN  Y4  2*  I )*/(?,«)*«(!  l4»(»,A»*ir(2l 
Y(NY43*  1 )*  7 ( I,  3>*W(  71471  l,»l*W(«) 

_ Y(  NY4A'  J > =M?  , ? >*W  ( ( > 47< 4 1 »W  (4  1 

4 5 N Y = N Y 4 J 

GF  T u r- M 
17  0 


^ 


Alternative,  but  less  efficient,  version  of  NEWVEC  which  better  Illustrates 
the  column  and  row  operations. 

SUBROUTINE  NEw'CCL(F,M  , N2 , 7 , nF^  N , A | ) 

DIMFNSION  VINM,NI,Z(7,(l,k(4l 

L_EERFQRM.__QLiK..K  _EXCHAM5E  JFANSECP  M U CN  Ch  CCLLMNS  Ml  T£  N CF  Y 


STARTING  WtTF  M|. 

*T*d  * NEWU  ,<t*IF*U  ♦ X*VI  ■ NEW  V 
N l , C 2 IS  N 2 BY  N 2 , AND  (l«V)"(Y(Ml)  , Y ( M 1 ♦ 1 I ,...)• 
C 2 ARF  STQRCC  IN  7.  F IS  SCALING  FACT  HR  FOR  X. 


! IN  THE  N14N2  ROWS  STARTING  W ITF  Ml. 

: COMPUTE  C2*(F*V  - *T*L>  * NEWU  ,<l*IF* 

: WHERE  Cl  IS  N 1 BY  N1,C2  IS  N2  BY  N 2 , AN 

: thf. matrices  a » c i , C2  arf  stqrcc  in  7. 

MsM I-  1 

DO  SO  K-Ml , N 
on  10  J= l ,N2 

V t J ).  - T lULTN-l  fcb&l  *E 

DO  10  L«l ,Nt 

10  WIJ ) *W( J )-Z(L ,J ) ♦ Y( M4L »K  I 

DO  20  J = 1 ,Nl 

U..M2 1 JO_  TIM t JL, KJL*F  ...  

DO  20  L*1.N2 

2 0 W(  N24  J)  = W ( N2  ♦ J ) ♦ 7 ( J,L)*Y(W4N1  ♦ L , K ) 

DO  J5  1 , N2 

S.«0^  _ ...  

DO  JO  L * 1 ,N2 


30 

35 

SaS+ZtJ,L4A)*W(L) 

Y<  M4 J ,K)  as 

D0  A5.Jal.Nl  . . . . 

S=0  . 

DO  AO  L a 1 ,N1 

AO 

SaS47(J,L+2)*W(N2+L) 

4 5 

Y(  M4N2F  J ,K1  *S 

50 

CONT INUE 

RETURN 

END 

SUBROUT INE  NEWROW(F,Nl,N2,7,Y,NM,N,Mt,t) 

DIMFNSICN  Y(NM,M  , 7(2,6  ),W(4  > 

PERFORM  BLOCK  EXCHANGE  TRANSFORMATION  ON  FIRST  1 ROWS  OF  Y 
IN  THE  NIFN?  COLUMNS  STARTING  WITH  COLUMN  Ml. 

COMPUTE  C2*IF*V  - XT*U>  ■ NFWU  ,C1*IF*U  ♦ X*V)  * NE  A V 

WHERE  Cl  IS  N1  RY  N1,C2  IS  N?  PY  N2,ANO  ( U , V ) = ( Y ( M l ),Y(M14I  ),  ...  ) . 

THE  MATRICES  X , C 1 , C 2 ARE  STORED  IN  Z.  F IS  SCALING  FACTOR  FOR  X. 

M = Ml  — 1 
DO  SO  K = 1 , I 
DO  10  J=1,N2 

JLLJ1  .......... 

DO  10  L= 1 , N 1 

[0  W<  J)»W(  Jl-ZIL,  J»*Y<K,MH.) 

DO  20  Ja l, NI 

_.  JX  ( N2EJ ) * YCK.MFJItF  

DO  20  L = 1 ,N2 

*0  W(N2HJ)=W(N2HJ)4Z(J,L)*V(K,MFN|FL> 

DO  35  J*1  ,N2 

•SaO. . 

00  30  L* 1 , N2 
JO  S«S4Z<  J,L4A)*WIL> 

JS  V(K,M4J)«S 

...  DO  . A S J*  1 , Ml . . ...  _. 

S»0. 

00  AO  La  I , N1 

10  S*S4Z<  J,L42»'">  AIN24LI 

IS Y IK_,J1*N2j4J  1*-S- ... 

>0  CONTINUE 
return 

END 
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